{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(200, 5695)\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAUMAAAD8CAYAAADt2MYTAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAIABJREFUeJztnXd0VVX6/p9NEhMhNKV3ghHsIIggFhRQsYuiCMsljooC6neJgqKOZcaxK/ZxxKHMiAqOUkQRsGDBAqEpozB0DCAl9IQW2L8/iP5y7vMCNzfJDTn3+azliufjTs4x771vzt377Pd13nsIIUSiU6GsL0AIIQ4HlAyFEAJKhkIIAUDJUAghACgZCiEEACVDIYQAoGQohBAAlAyFEAJAMZOhc+5C59xC59xi59x9JXVRomxRXMOLYntgXKw7UJxzSQD+B6ALgGwAMwFc573/+UDfk5SU5JOTkwPu6KOPpnF5eXnkKlasSC7yZx1o3L59+8itXbuWXOXKlQ95Lbm5udi1a5ejgSEhlrg657xzwV9J3bp1ady2bdvIVapUiVxKSgq5tLQ0ctZrd+XKleRq1apFznqNbdy4cYP3vib9h5BQ1NgmJyf71NTUgKtXrx6Ny8nJIVejRg1ykT8LAPbs2UNu79695FavXk2uatWqUX3vunXrooorZ5PoaQtgsfd+KQA4594FcDmAA75pkpOT6ZfZu3dvGjdnzhxyLVu2JFe9enVyrVu3Jme98F944QVyHTt2JDd79uzA8ZQpU2hMyChyXJ1z9ELv27cvjfv888/JtWvXjpyVSFu0aEFu165d5O644w5yt99+OznrNTZq1KgVJMNFkWKbmppKv/fHHnuMxo0YMYLczTffTC4jI4OcleSsP5p//vOfyV100UXktm7dSu6ll16KKq7F+ZhcH8CvhY6zC5wo3yiu4UWxPQjFuTO0PirS5xbnXB8AfQAgKSmpGKcTcaLIcRXlhkPGtnBcjzjiiHhc02FDce4MswE0LHTcAADd83rv3/Det/Het1EyLBcUOa6R84XisOWQsS0cV2tOPswU5/92JoBM51xTAKsA9ADQ82DfkJKSgvr1g3flmzdvpnETJ04kt2IFf+zv2rUruUWLFpGbPHkyOWtusXnz5uSWLVsWOE6AF0iR41qxYkWccsopAbdhwwYa9+2335Jbs2YNubPPPptco0aNyI0cOZJcnz58sxp5bQCQlZVFLgEoUmxTU1PRrFmzgFuyZAmNs96b48ePJ5eenk6uTZs25KZPn07u3nvvJbdp0yZy77zzDrloifmd7b3Pd87dDmAygCQAw7z3/435SsRhgeIaXhTbg1Os2xzv/ccAPi6haxGHCYpreFFsD4x2oAghBJQMhRACQDE/JheVvXv30oLJCSecQONGjRpFznrY89xzzyVnTdJbD2dnZ2eT++STT8g1aNAgcGw9MS/2P3hdmMiJdwAYM2YMuf79+5Pr3r07uVmzZpFr1aoVucgFL8CO63HHHUdOBMnLy6NNB5HvBwB46623yB1//PHkrEWV999/n5y1e8XKCdaOltNOO43cl19+Sc5Cd4ZCCAElQyGEAKBkKIQQAJQMhRACQDFKeMVC48aNfeST5CeeeCKNmzZtGrmGDRuSmzFjBrlrrrmG3OOPP07O2pVgVdX49NNPA8dffPEFNm3apP1nhcjMzPQvvvhiwFm7Df77X36+16pQ8ssvv5Czqhs9+uij5KwSbldddRW5r7/+mtzw4cNnee95S0SCkpmZ6Z9//vmAW7duHY2zYmhVozn11FPJWa8Ta7H0jDPOIGe9X999911y06dPjyquujMUQggoGQohBAAlQyGEAKBkKIQQAOK8A2Xr1q20IPHee+/RuIEDB5J79tlnoxrXr18/cpdeeik5q9/JN998Qy6yZLy1CJDorF+/Hv/4xz8Cziq/fuedd5IbO3YsuQcffJDcOeecQ85aQIl2oa1Xr17khg8fTi6RseJqFXy1FkbmzZtHrmdPrhZ22WWXkevcuTO5Jk2akLNK80W7U81Cd4ZCCAElQyGEAKBkKIQQAIo5Z+icWw5gG4C9APL1wGp4UGzDieJ6YEpiAeVc7z03vDDYt28fdu7cGXCdOnWicc888wyfxJgEbdq0KblBgwaRsxZVrIl7qwH5K6+8Eji2nsAPMVHFNjk5GUcffXTAWT1mvv/+e3LXXnstud9++43cxx9zcebBgweTGzJkCLnXXnuNXIITVVwrVKhAO0TOPPNMGvfVV1+Rs3oaW/G/5557yI0bN46cdV6rifwHH3xALlr0MVkIIVD8ZOgBTHHOzSrotyrCg2IbThTXA1Dcj8kdvPernXO1AEx1zi3w3gfumQs3pU5LSyvm6UQcOWhsC8e1UqVKZXWNouhEHVer6EWYKdadofd+dcHXdQDGAmhrjPmjKbX1wKY4PDlUbAvHVX/kyg9FiatVVj/MxHxn6JyrBKCC935bwb+fD+AvB/uezMxMmgi3epZYTaTXr19Pzppo/fzzz8lZZX0mTZpE7oknniAXWYpo4cKFNCZsFDW2TZo0wbBhwwLOavBu/e4mTpxIrl27duRGjBhB7tVXXyX32WefkRs6dCi57du3k3v99dfJhYmixjUjI4PeO1a/k127dpFbvHgxOathvLXr680334zq53Xo0IGctdvs5ZdfJmdRnI/JtQGMLWgElAzgbe89d94R5RHFNpworgch5mTovV8KgCukinKPYhtOFNeDo0drhBACSoZCCAEgziW8li5dih49egTc3LlzaVydOnXIWRPto0ePJmeV2IrszwEAV199NbkLLriAXGQJr9zcXBqT6Cxfvhw33HBDwH3xxRc0rl69euQyMzPJWf0urJJQf/7zn8n93//9X1TjbrnlFnIiyJo1a6gk1ief8BRjzZo1yf3nP/8hZ/WiqVWrFjmrDNd1111HbsKECeR27NhBLlp0ZyiEEFAyFEIIAEqGQggBQMlQCCEAxHkBJS0tDc2bNw+4yKbygF0SaMWKFeT+9a9/kdu0aRO5tm1plyD27NlD7q677iK3du3awPHevXtpTKKTnJyM2rVrB5y1Y8TaWWQ1jLcWxqxG5dYi2Jo1a8hZJdysMnEiSH5+Pu38st4j1iKo1bMmshweADRq1IhccjKnJasp/VlnnUWufv365KJFd4ZCCAElQyGEAKBkKIQQAJQMhRACQJwXUHbt2oUlS5YEXEEFjQBJSUnkPvzwQ3J33303OatRuVXq6f777ydnlQ7Ky8sLHO/bt4/GCI6jVbvSiutPP/1EbsCAAeSsXinvvPMOuTlz5pDr3r07OavxubVrJtGpUCF4v9S+fXsaYy2CWWXzrrjiCnJWbIYPH07OWrh54YUXyFmvsWjRnaEQQkDJUAghACgZCiEEAMB57w8+wLlhAC4BsM57f2KBOwrAaABNACwHcI33np92jqBBgwY+ck5vwwZu32r1L7YqnsyePZvc8ccfT87quWr17bAq40Q+7Pndd99hy5YtPNFZDimp2NavX9/3798/4KyHZK3f+ebNm8mtWrWKXGRfZgA48cQTyWVkZJCbNWtWVO6zzz6bFYam6iUV19q1a/vIajHWnPmWLVvIWXG1Ksqcd9555KxNEtaD2KNGjSJntQcZO3ZsVHGN5s5wBIALI9x9AD7z3mcC+KzgWJQ/RkCxDSMjoLgWmUMmw4I2ghsj9OUAfu/4MxIALxOJwx7FNpworrER65xhbe/9GgAo+MoVGgtwzvVxzmU557JUGLVcEFVsFddyR5HjWpxCqeWRUl9AKdyHVc3Gw4PiGk4Kx/XII48s68uJK7E+dL3WOVfXe7/GOVcXwLpovmnbtm2YNm1awFlvJOth6rFjx5Jr0KABOevBTqucv1XxxGo3EFktJQGq1hQ5tmvXrsXzzz8fcN26daNx11xzDbkZM2aQ27lzJ7mbb76Z3LnnnkuuRYsW5KyFluI8nFtOKXJc8/LyaKHp9NNPp3FWv/Hx48eTs1pydOnShZzVX9lq52B9IrHaQ0RLrHeGEwD83vTiBgD8fy7KK4ptOFFcD8Ehk6Fz7h0A3wFo7pzLds7dBOBJAF2cc4sAdCk4FuUMxTacKK6xcciPyd57bku1n04lfC0izii24URxjQ3tQBFCCESxA6UkqVatmj/77LMDrnr16jTOuibLffbZZ+SeffZZctYk/cKFC8lZZcRPOumkwPGAAQOwaNGiUOxAKSlq1qzpIxeucnJyaJy1+GSVaf/xxx/JPfroo+Sstg/r1vG6QNeuXclFtp8AgAsvvDAUO1BKiho1avjLLrss4Hbv3k3jqlSpQi4lJYWc1c5jzJgx5Pr06UPOag9w+eWXk7P6MJ911lkltgNFCCFCj5KhEEJAyVAIIQAoGQohBIAy6Jt87LHHBpxV9v+qq64it3z5cnJWuZ7IHS4AcOutt5Kz+utak/mRZYKs/r2JjnOOSiydfPLJNM6aGLcWS7Zv307OKvE+ZMgQcla/ZmuXg1X+SQRJT0+nRUXrfWOV4bL2Nf/888/kBg4cSM5qGWC1ZGjcuDE5K/7RojtDIYSAkqEQQgBQMhRCCABKhkIIASDOCyg7duygMj7WooW128TqzXrcccdFdd4ePXqQO+qoo8hZpYi+/PLLwLEWUJjc3FzMnDkz4Ky+xPn5+eQmTJhAztpFkJ6eTq5z587katSoQc4qHbZ06VJyIsj27dvp9b9r1y4ad/3115N7+umnyV155ZXkrDJcHTp0IBe58ArYJdysXtrRojtDIYSAkqEQQgBQMhRCCADRFXcd5pxb55ybX8g94pxb5ZybW/DPRaV7maKkUVzDi2IbG9EsoIwA8AqAyHpJQ7z3XC/rIKSmpqJZs2YB17JlSxpn7Syw+p1Yixm9evUil5mZSc4qO/Tbb7+Ri3y6/tdff6Ux5ZQRKKG4Vq1aFZdccknAWYtbt912Gzmrp8Z3331H7sYbbyRn9aw55phjyFkN46tWrUouRIxACcQ2LS0NJ5xwQsBZi1uDBg0iZy2WTJ8+ndw999wT1bVYsbYW1R544AFyAwYMiOocsfZNFuUcxTW8KLaxUZw5w9udcz8W3JJzhVZRXlFcw4tiexBiTYZ/B9AMQEsAawA8d6CBidyUuhwSU1zz8vLidX0idqKKbeG4Ws8AhpmYkqH3fq33fq/3fh+AoQDaHmRswjalLm/EGteKFSvG7yJFTEQb28JxtXqah5mYdqD83oy64PBKAPMPNv53rJ0KVm+DKVOmkLOeaLf6XaxcuZLcsGHDyK1du5bc0KFDyc2ePTtwHOa7oFjjmpeXh6ysrICzmrnPn88/LrInDgBccMEF5FavXk1u6tSp5KwdLffddx+5ZcuWkQszscR2586dVHYrckEFAL2nAXthtF69euSs95O1gGa9Xx966CFyr7/+OrloOWQyLOjB2hFADedcNoCHAXR0zrUE4AEsB8AFA8VhjeIaXhTb2Ii1b/I/S+FaRBxRXMOLYhsb2oEihBBQMhRCCABxLuFVrVo1XHzxxQH32muv0bhbbrmFnNUDxSrD9corr5CzSj1ZE7cff/wxuXvvvTdw3K9fPxqT6Fi9Mnr37k3j+vbtS27u3LnkrF0p1mKJ1QOjQgX++x5ZNg4A+vfvT27y5MnkEpn09HQqp2XF8JFHHiFn7fo54ogjyFnl9aweSNbC6NatW8mNHDmSXPv27clZ6M5QCCGgZCiEEACUDIUQAoCSoRBCAIjzAsqWLVswadKkgLNKKVmN5atVq0bulFNOIZednU3OWiyxzmG50aNHB443bdpEYxKdnJwcat69d+9eGle9OtcGsMqrWeW/IncCHQjrdWKdw5poF0FycnLw1ltvBZy1aGntYbZKqVnO2pVibQO04rpw4UJyxYmr7gyFEAJKhkIIAUDJUAghACgZCiEEgDgvoFg7FawyXNaT6lbz6rfffpvcDTfcQM7q23DGGWeQsyZpI3tvpKSk0JhEp0qVKjj//PMDrlWrVjQuIyOD3Lx588gtXryYXM+ePcnVr1+fXGQvlgP9vE6dOpF7//33ySUyycnJtGBi7fqZMWMGOSsO1m4jayfQDz/8QO6ii7h/VceOHcn96U9/IvfRRx+Rs9CdoRBCQMlQCCEAKBkKIQQAJUMhhAAAOO99/E7m3HoAKwDUALChlE9XWudo7L2vWQo/t9xSKK5A6ce2NH++YluIRItrXJPhHyd1Lst736a8n0Mwpf17V1zLhkSIqz4mCyEElAyFEAJA2SXDN0JyDsGU9u9dcS0bQh/XMpkzFEKIww19TBZCCCgZCiEEgDJIhs65C51zC51zi51z95XSOZY7535yzs11zmWVxjlEEMU1nCRSXOP90HUSgP8B6AIgG8BMANd5738u4fMsB9DGe1/aD3YLKK5hJdHiGu87w7YAFnvvl3rvdwN4F8Dlcb4GUfIoruEkoeIa72RYH8CvhY6zC1xJ4wFMcc7Ncs71KYWfL4IoruEkoeIa1+KuALj93P5fREnTwXu/2jlXC8BU59wC7/1XpXAesR/FNZwkVFzjfWeYDaBhoeMGAFaX9Em896sLvq4DMBb7b/dF6aG4hpOEimu8k+FMAJnOuabOuSMA9AAwoSRP4Jyr5Jyr/Pu/AzgfwPySPIcgFNdwklBxjevHZO99vnPudgCTASQBGOa9/28Jn6Y2gLEFDeGTAbztvf+khM8hCqG4hpNEi6u24wkhBLQDRQghACgZCiEEgGImw3hs1RHxR3ENL4rtgYl5zjCWrTpJSUk+OTm4ZtOwYUMal5OTQ65u3bpRXZfVbL5CBc75y5YtI3f00Ucf8ufl5eVh9+7d1vNXoSCWuFaoUMEnJSUFXL169Wjc9u3byUU2KQeA3bt3kzviiCOs85LLzs4mV7VqVXJ79+4lt27dug1h7oFS1NimpKT41NTUgLPiunnzZnLVqlUjt2/fPnKRP/9A49avX08uLS2NnMWqVauiimtxVpP/2KoDAM6537fqHPBNk5ycTEntueeeo3EjR44k99BDD0V1UYsWLSKXnp5OrlevXuSuvfZacitWrAgcf/nll1FdRzmmyHFNSkqipDZ48GAa980335Dr0aMHudWr+VE264/mkUceSW7gwIHkLr74YnJbt24lN2TIkBUkw0WRYpuamopTTjkl4B599FEa9/7775Pr1q0buR07dpDLyMiIatxrr71Grnnz5uQsBg8eHFVci/MxOaqtOs65Ps65LOdclvXXWBx2FDmu1l9ycVhyyNgWjuuePXvienFlTXGSYVRbdbz3b3jv23jv20R+lBKHJUWOq/VxVRyWHDK2heOakpISp8s6PCjOqzguW3VE3FFcw4tiexCKM2f4x1YdAKuwf6tOz4N9Q2pqKjIzMwMuK4trOVqTpdYcVLt27chZC0J5eXnkBg0aRG7nzp3k5s8P7gxKgI+ERY5reno6zjrrrICzJtWnTZtGrmLFiuSseaQaNWqQ++mnn8hdfjlXmLJeJ88//zy5BKBIsU1NTUXjxo0DbvLkyTTugw8+ILdlyxZyXbt2JWfN3X788cfkrDnjVq1akRs1ahS5aIk5GcZpq46IM4preFFsD06x9iZ77z8GwGlclGsU1/Ci2B4YzXwLIQSUDIUQAkCcS3jt2rULS5cuDbhLLrmExlkPWLZu3Zqc9WCnNcF78sknkxs/fjy5KlWqkIucQP7xxx9pTKKzc+dOLFiwIOCsRYtx48aRO/PMM8kNHz6c3M8/83PB+fn55CKvAwC+/vprctbradKkSeQSmR07dtDvvXv37jTOctYmiRYtWpD79ddfydWpU4ectbPoxRdfjOoc0aI7QyGEgJKhEEIAUDIUQggASoZCCAEgzmX/W7Ro4d98882A+/bbb2lcbm4uOWt3yKmnnkouskQYADz22GPkmjZtSq5Ro0bkIiut/PLLL8jNzQ1tCa9YyMjI8I8//njAWXGwyqZZiyBWuS4rNsOGDSNnTb5bk+rWhPwrr7wyy3vfhv5DgtKiRQs/dOjQgLPKof3www/kGjRoQG7JkiXkOnToQC7ynABQu3ZtctZrYtasWeReffXVqOKqO0MhhICSoRBCAFAyFEIIAEqGQggBIM4LKOnp6f6kk04KuJYtW9K40047jdwnn3Bf6b59+5I777zzyD311FPkrN4bVomp2267LXB88803Y8GCBVpAKUTNmjV9ZOksaxGsX79+5KzdJnfccQe5008/ndzTTz9NzupjM2bMGHIPPPAAufbt22sBpRA1a9b0V1555SHH3Xcf95WKXFAD9r93IrFibbWCiGw/AABTpkyJ6hzNmzfXAooQQkSLkqEQQqCYhRqcc8sBbAOwF0C+PmKEB8U2nCiuB6Ykqtac673fUAI/Rxx+KLbhRHE1iGsJr7S0NBx33HEBZz1ZbvVAaNOG/4BZ/VNefvllclYf5rvuuouctRvib3/7W+B4zZo1NCbR2bdvH3bt2hVwVryshaxjjz2W3NixY8lZPUusuN5+++3krP4p1sKNCFKhQgXqTW3t3Lr//vvJde7cmVz16tXJ/eUvfyFnLb5YvW02bOB8/u6775KLluLOGXoAU5xzs5xzfYr5s8ThhWIbThTXA1DcO8MO3vvVzrlaAKY65xZ4778qPKDgF94HACpVqlTM04k4ctDYKq7llqjjWrly5bK6xjKhWHeG3vvVBV/XARgLoK0x5o+m1GlpacU5nYgjh4pt4bimpqaWxSWKGChKXCM/IoedmJOhc66Sc67y7/8O4HwA8w/+XaI8oNiGE8X14MS8A8U5l4H9f1mA/R+33/be/+0g34I2bdr4yEWPDz/8kMZZvSispuTnnHMOudGjR5OzeiXMmTOHXN26dcl16dIlcNymTRtkZWWFegdKUWPbunVrH1nGad68eTTOmizftm0bOasvzksvvUTO2jFkndfavWDd9dSpUyfUO1CKGlfr/Wr1ALJKbu3evZucVXLP6llk7Q6aPXs2uUsvvZRcXl4euWbNmkUV1+I0kV8KgF9lotyj2IYTxfXgaAeKEEJAyVAIIQAoGQohBIA470BZtWoVHnzwwYCzdiXccMMN5KwJ1FtvvZXc8ccfT27IkCHkIhdGALucUOTT9Tk5OTQm0VmyZAmuuuqqgLOavrdq1YqctYvo73//Ozlrt8Ett9xCrlevXuSsMlQ33XQTORFkxYoVVMJu5syZNO6yyy4jN27cOHLWAprVF2XUqFHkrEWwG2+8kVzv3r3JRYvuDIUQAkqGQggBQMlQCCEAKBkKIQSAOPdAqV27tu/Zs2fAnX/++TRu+vTp5GrWrEnOenrdahhu7V6xmmFbuxLq1asXOH766aexcuXKUO9AKSq1a9f2kQsX7dq1o3FLly4ll5SURG7FihXkrBgec8wx5Bo3bkxu7ty55Dp27EiuW7duod6BUlSaNGniIxc8TzjhBBq3YMECcuvWrSNnLYJVrFiRnPVe37dvH7mMjIyoznvTTTepB4oQQkSLkqEQQkDJUAghACgZCiEEgDjvQMnPz8f69esDzlpA+f7778l9+umn5M444wxyF154Ibnx48eTs3olWA2oIxdarAn/RCcpKQlVqlQJuNatW9O47777jpy1A+WKK64gd9JJJ5GzYvjss8+Se+KJJ8jt2bOHnAiyc+dOWhyxdptYu0isWFvvr8ieSAAwZsyYqNxf//pXcs2aNSMXLbozFEIIKBkKIQSAKJKhc26Yc26dc25+IXeUc26qc25RwVfuASgOexTbcKK4xsYhH7p2zp0NYDuAf3nvTyxwTwPY6L1/0jl3H4Dq3vt7D3Wyo446ykf2U7WaCf3222/kWrZsSc7qr2z1V7UqXlhzC1Yv5cj5kA0bNmD37t2heOi6pGJbs2ZNb1WGiWTTpk3krJ7GVmn59u3bkzv99NPJ1alTh9y///1vcsuWLSP36aefhuKh65KKa5MmTfzDDz8ccD/99BON27lzJ7latWqRsx7Otlp3WO916+Fsq7rNL7/8Qm7ixIkl89B1QRvBjRH6cgC/d/AeCYBnvMVhj2IbThTX2Ih1zrC2934NABR85T8Doryi2IYTxfUQlPoCinOuj3MuyzmXtWvXrtI+nYgTheNqfUwS5ZPCcd2+fXtZX05ciTUZrnXO1QWAgq+8O7oANRsvd0QV28JxTUtLi+sFipgoclzT09PjeoFlTawPXU8AcAOAJwu+8lPNBvn5+TSJ3qFDBxp39913k5s6dSo5a1Ldeoi7U6dO5Pr06UOucuXK5CpVqhQ43rgxciomdBQ5tnl5edSWwaoe9OSTT5Kz+iFbD2xbcb3gggvInX322eQiKw8BwKxZs8iFnCLHdcuWLdTX3FqMvOaaa8hZlaesijdnnXUWOesB+3vuuYfc3r17o3LREs2jNe8A+A5Ac+dctnPuJuz/hXZxzi0C0KXgWJQzFNtworjGxiHvDL331x3gP/HtlihXKLbhRHGNDe1AEUIIKBkKIQSAOFetcc5R1Zc5c+bQuPnz55PbvXs3uTVr1pC7+uqrybVt25actRvC6pucmZkZOI5mp0WikZKSgvr16wecc7xJZ/DgweSsJwy++OILchdddBE5a0J+0aJF5F544QVyDzzwADlr4j6RSU5OxlFHHRVwixcvpnHbtm0jZ71fZ8yYQc6qgmO9TlatWkWuW7du5Kz356RJk8hZ6M5QCCGgZCiEEACUDIUQAoCSoRBCAIjzAkp6ejrtEFi7di2N6927Nzmrl+5TTz1F7o033iA3fPhwch999BG5hg0bkossE5STk0NjEp3KlSvj3HPPDTirf621CLJ8+XJyVqknq8S7tdA2bty4g13qH7z99ttRjUtkKlWqRGXSrIUMq5x/dnY2OWuvs1VezVrIXLlypXl9kXz++efkokV3hkIIASVDIYQAoGQohBAAlAyFEAJAnBdQcnNzqZ9qSkoKjatZsya5YcOGkevevTu51atXRzUuIyODnLWY8+KLLwaO1W+XycnJwciRIwOuSZMmNG7AgAHkrD42Vgkva0fLMcccQ65atWrkrruO6xa8+uqr5ESQrVu30oJEo0aNaNzmzZvJWSW8rPecVRKvf//+5E488URyHTt2JDdlyhRy0aI7QyGEgJKhEEIAUDIUQggAsTdGt/FJAAAF0ElEQVSRf8Q5t8o5N7fgH36aVhzWKK7hRbGNjVibyD8CYLv3/tminKxBgwb+zjvvDDir2bRV6sfqd2GV/7L6MUybNo2ctdskPz+fXOST9M888wxWrlxZ7pvIl2RcGzVq5CP71lg7EKyyaVbDcKuLotXHxCrhZi2gRJZhA+wyYQ8//HCYm8g/giLGtlmzZv7xxx8POGsHilWuKzmZ12atxcczzzyTnBXX2rVrR3XeDRs2kOvZs2epNpEX5RzFNbwotrFRnDnD251zPxbcklcvsSsSZY3iGl4U24MQazL8O4BmAFoCWAPguQMNLNyUOjc3N8bTiTgRU1wTrdl4OSWq2BaO69atW+N5fWVOTMnQe7/We7/Xe78PwFAAXFf//4/9oym1VWVCHD7EGtdEazZeHok2toXjWqVKlfheZBkT0w4U51xd7/3vs5xXAuBaSgbbtm2jJ9p79OhB4yZOnEjOmhi3JkutJ9qzsrLIDR06lNxDDz1ELrLJuVWaKizEGteNGzfinXfeCThr14e1i6h6df60ZrkaNWqQi9wdBAALFy4kN2HCBHL//Oc/yYWZWGK7ceNGjB49OuCuv/56GmeVzbN6m1h/NK0F1Pfee4/ct99+S27EiBHkIkvuFYVDJsOChtQdAdRwzmUDeBhAR+dcSwAewHIAt8Z8BaJMUFzDi2IbG7E2kU+sP6shRHENL4ptbGgHihBCQMlQCCEAxLmEV6VKldC+ffuAu+mmm2jcc8/xqr+1K8Uq9WT9vK5du5KzShFZT8i//vrrgeN+/frRmESnevXquPbaawPOWoyySqn9+uuv5KpWrUqub9++5Nq25QXRFi1akLMan1s7UJo1a0YukalatSouvvjigBs0aBCNs0qzWY/RJSUlkevTpw+5Sy+9lFydOnXIWWXCbrvtNnKR5eUOhO4MhRACSoZCCAFAyVAIIQAoGQohBIA4L6Ds2LED8+bNC7h27drRuMaNG5OztgZZT68PHDiQnLUrJS0tjZzVqPqtt94KHFt9UhKd3NxcWuCy+thYcbV2lhx55JHkrB0tVlkna0Fm2bJl5KwFFBHEer9a8apYsSI5K/7169cnZy2WWO9Nq7eRtSAzdepUctGiO0MhhICSoRBCAFAyFEIIAEqGQggBIM4LKPn5+Vi/fn3AnXzyyTRu+PDh5K688kpy1hPo1gLK119/TS5yxwQAnHPOOeQie69YfVcSndTUVDRt2jTgtmzZQuOs3iYVKvDf4x9//JHcgw8+SM5aGLN2pUReGwDaWSGYlJQUWvSweptYZdOsxTJrcatz587kmjRpQq5Dhw7kjj32WHKdOnUiFy26MxRCCCgZCiEEACVDIYQAoGQohBAAomgiX6Inc249gBUAagDgBiYlS2mdo7H3nh+vT2AKxRUo/diW5s9XbAuRaHGNazL846TOZUXT4f5wP4dgSvv3rriWDYkQV31MFkIIKBkKIQSAskuG3Gi1fJ5DMKX9e1dcy4bQx7VM5gyFEOJwQx+ThRACZZAMnXMXOucWOucWO+fuK6VzLHfO/eScm+ucyyqNc4ggims4SaS4xvs5wyQA/wPQBUA2gJkArvPe/1zC51kOoI33vrSfZRRQXMNKosU13neGbQEs9t4v9d7vBvAugMvjfA2i5FFcw0lCxTXeybA+gMJ1fLILXEnjAUxxzs1yznGXalHSKK7hJKHiGtd6hgCc4Urjc3oH7/1q51wtAFOdcwu891+VwnnEfhTXcJJQcY33nWE2gIaFjhsA4LZXxcR7v7rg6zoAY7H/dl+UHoprOEmouMY7Gc4EkOmca+qcOwJADwATSvIEzrlKzrnKv/87gPMBzC/JcwhCcQ0nCRXXuH5M9t7nO+duBzAZQBKAYd77/5bwaWoDGOucA/b//73tvf+khM8hCqG4hpNEi6t2oAghBLQDRQghACgZCiEEACVDIYQAoGQohBAAlAyFEAKAkqEQQgBQMhRCCABKhkIIAQD4f9Vxs2HTYy1iAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x21d00023ef0>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "from skimage import io\n",
    "import matplotlib.pyplot as plt\n",
    "import numpy as np\n",
    "import os\n",
    "import pandas as pd\n",
    "\n",
    "def randomfaces(d,m):\n",
    "    ## create a matrix d*m to hold the resault\n",
    "    matrix_one = np.empty([d,m])\n",
    "    ## the random matrix \n",
    "    matrix = np.random.randn(d,m)\n",
    "    ## create an empty array\n",
    "    l2 = np.empty([d,1])\n",
    "    ## compute the l2 norm in each row\n",
    "    for i in range(np.shape(matrix)[0]):\n",
    "        l2[i] = np.linalg.norm(matrix[i,:],ord=2)\n",
    "    ## divide each row by l2\n",
    "    for i in range(np.shape(matrix)[0]):\n",
    "        matrix_one[i,:] = (matrix[i,:]/l2[i]).reshape(-1)\n",
    "    ## aha!\n",
    "    return matrix_one\n",
    "\n",
    "def showEigFace(Eig, row, col):\n",
    "    plt.figure()\n",
    "    for index in range(row * col):\n",
    "        plt.subplot(row, col, index + 1)\n",
    "        plt.imshow(Eig[:, index].reshape(10,20).T, cmap=\"gray\")\n",
    "    plt.show()\n",
    "    \n",
    "if __name__ == \"__main__\":\n",
    "#     faces_cls = os.listdir(root)\n",
    "#     faces_img_path = [os.path.join(root, cls, '1.pgm') for cls in faces_cls]\n",
    "    \n",
    "#     origin = io.imread(faces_img_path[0])\n",
    "    \n",
    "#     print(origin.shape)\n",
    "    \n",
    "    test = pd.read_csv('test.csv')\n",
    "    test = np.array(test)\n",
    "    test_ones = np.ones((5695,), dtype=int)\n",
    "    test_new = np.row_stack((test, test_ones)) \n",
    "    test_new = np.array(test_new)\n",
    "    \n",
    "#     print(origin.shape[0]*origin.shape[1])\n",
    "#     X = make_input(faces_img_path)\n",
    "#     print(X.shape)\n",
    "#     matrix_one = []\n",
    "#     matrix_one = np.array(matrix_one)\n",
    "    matrix_one = randomfaces(200,1024)\n",
    "    eig_vec = np.dot(matrix_one,test_new)\n",
    "    print(eig_vec.shape)\n",
    "    showEigFace(eig_vec, 2, 3)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
